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Abstract 



We study the coefficients of the expansion F{R) = ^ c^R^ + \ C2-R^ + ciR of the free 
energy of spherical bubbles at T = Tc in pure glue QCD using lattice Monte Carlo 
techniques. The coefficient C3 vanishes at T = Tc and our results suggest that the sign 
and the order of magnitude of ci is in agreement with the value ci = ±327rTc/9 (- for 
hadronic bubbles in quark phase, + for quark bubbles in hadronic phase) computed by 
Mardor and Svetitsky from the MIT bag model. The parameter C2 is small in agreement 
with earlier determinations. 



1) University of Helsinki, Department of Theoretical Physics, Siltavuorenpenger 20 C, 
00170 Helsinki, Finland; kajantie@finuhcb. 

2) Research Institute for Theoretical Physics, Siltavuorenpenger 20 C, 00170 Helsinki, 
Finland; permanent address Universitat Bielefeld, Fakultat fiir Physik, Postfach 8640, 
D-4800 Bielefeld 1, BRD; leo@c240.uni-bielefeld.de. 



3) CERN/TH, CH-1211 Geneve, Switzerland; rummukai@cernvm. 



- 1- 



Bulk QCD matter is stable in the quark-gluon plasma phase (Q) for T > Tg and in the 
hadronic phase (H) for T < Tc. At T = Tc the free energies coincide and there is a 
latent heat L: 

F,{T,) = -Vp,{T,) = Fh{T,) = -Vphin), L = T,\p'^{T,) - p'f^{T,)], (1) 

assuming the transition is of first order and taking the chemical potential n = 0. The 
determination of the free energy has been studied in great detail with lattice Monte 
Carlo techniques [1] . At T = Tg a stable H-Q (order-disorder) interface can exist in the 
system and contributes to the free energy: 

F = -pV + aA. (2) 

Lattice Monte Carlo computations [2-3] using planar interfaces have given the value 

a(Te) ^ O.IT^ (3) 

for the interface tension in the case of pure glue. 

The purpose of this letter is to go still one step further in the lattice Monte Carlo 
determination of free energies of domains of one phase immersed in the other: spherical 
bubbles of H phase in Q matter (or Q bubbles in H matter) [4] . Then the free energy of 
an H bubble (relative to a system entirely in Q phase) with radius R can be expanded 
in powers of as follows: 

F{T, R) = [p,(T) - phiT)] IttR"" + a(T)4rf + ^{T)R +... (4) 

The new element here is the last curvature term. Using the MIT bag model, i.e., treating 
the bubble as a spherical cavity, gluons as eight copies of Abelian photons, and assuming 
bag boundary conditions, Mardor and Svetitsky [4] found that 

«bag(T) = 0, ^bag(T) = -- ^T^; (5) 

for Q bubbles in H the sign of 7 is positive. Our main result is that also lattice Monte 
Carlo shows evidence of the pattern in eq. (5) : a small area term a and a large curvature 
term 7, the numerical value of which is compatible with that in eq. (5). 
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If only the terms shown in (4) contributed to F{R), there would be an -Rmin with 
F{Rmin) < F{0) = 0, i.e., H bubbles would be stable even for T > Tg. However, even 
in the bag model the terms in 1/R omitted in eq. (4) decrease F{0) — F{Rmin) and on 
the basis of our results presented here one cannot conclude anything on the stability of 
H bubbles for T > Tc. 

Before presenting the lattice results, consider eq. (4) phenomenologically for T Tc. 
Scaling with T^: f = T/T^, R = T^R, we firstly can write Pq{T) - ph{T) = L{f - 1) and 
estimate the latent heat by L = TcPq{Tc) = 327r^T^/45. For H-Q interfaces we do not 
know a{T) away from Tg, but a study of order-order interfaces in the Q phase [5] has 
given a very large value, Ta'{T)/a{T) 30, for the surface entropy, Sg = —a'{T)A. 
We use this same value and write a{T)/T^ = 0.1 + 3(T — 1). Our lattice data below 
indicates that 7 depends only weakly on T and we finally have 



Quantitatively, this model implies that at T = « 200 MeV a stable H bubble has a 
radius of about 4/Tc=4 fm. The radius has shrunk to 1 fm at T = l.OSTg. 

Our aim is to determine the coefficients of the expansion (4) at T = Tg. First, we 
point out that the expected magnitude of 7 = 7/T'^ = — 327r/9 is much larger than 
a — a/T^ = Q.l. This diff'erence is essential for the success of our calculations: as it will 
turn out, the actual thickness of the bubble wall is of order l/T^, so that the definition of 
the radius of the bubble is necessarily ambiguous. For example, possible non-equivalent 
choices are the radius where the order parameter reaches some fixed value between the 
ordered and disordered states; the turning point of the radial order parameter; or the 
corresponding quantities for the energy density. 

This ambiguity affects the functional form of F{R). Let us redeffne the radius R 
R^ = R- 5, so that at T = T^. 



F{T,R)/Tc = 47T 




(6) 



F = A-Kaff + 7^ = ATTaRl + (7 + 87Ta5)Ri + 



(7) 



If now a and 7 are of the same magnitude, even a slight variation in S will spoil the 
0{Ri) term. In fact, in this case a prime candidate for the bubble radius would be a 
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value of R where F = AiraR^ exactly! However, the values given above for 7 and a 
imply that 5 would have to be Ki A/T for the linear term to vanish. This is more than 
the wall thickness W ~ 2/T, so that this 5 would move the radius completely outside of 
the actual interface. Thus, the whole concept of the curvature contribution to the free 
energy makes sense only when 7 S> SiraW ^ Sira/T. 

We would like to stress that this ambiguity in the subleading terms is of quite general 
nature; it is present both in physical phenomena and in lattice calculations. We shall 
return to this question more in detail further down. 

On a homogeneous system at thermal equilibrium the bubbles are rare and lattice 
calculations would be extremely time-consuming. To remedy this and to gain systematic 
control of the bubbles, we choose a spatially spherical region of radius R and adjust 
the couplings /3 = 2Nc/g'^ = /3c ± A/? in the interior and exterior, respectively. This 
acts like an external field and stabilises the bubble. For simplicity, we call this R the 
radius of the bubble, and check the sensitivity of the results to reparametrisations later. 
The coefficients of eq. (4) are now evaluated with various A/3, and the final answer is 
obtained by extrapolating A/3 0. In the calculations one cannot use too small a 
A/3, otherwise the bubble will vanish. We assume that the relevant quantities behave 
linearly in A/3, as long as the bubble is preserved. 

It is not possibly to evaluate directly the free energy F using a MC simulation. However, 
the derivatives of the free energy are accessible as ensemble averages, thus providing, 
after integration, the value of the free energy up to a constant. The 'standard' way to 
do this is to integrate with respect to the inverse coupling /3 [1,2,3]. However, as we are 
interested in the i?-dependence of the free energy, it is natural to take the derivativee 
w.r.t. R, the bubble radius. With this method we avoid the very delicate and time- 
consuming integration over /3. However, in order to measure the derivative, we have to 
modify the action to form 

V/3(rn, R)S{n) = V(/3e + A/3tanh !i^)5wiison(n), (8) 

z — / z — / uu 

n n 

where the sum is over the lattice, Vn is the distance from the center of the bubble and w 
(chosen as 0.5 lattice spacings a) parametrises the width of the cross-over of the function 
/3(r, R); w is not to be mixed with the actual thickness of the interface. The quantity 
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accessible by means of MC then is 



d{F/T) 



dR 



A/3 



w 



/^cosh-(!°-^)S(„)\. 



(9) 



where the average is taken over the ensemble created with the action of eq. (8). 

The role of the parameter w deserves some elaboration. It was introduced to facilitate 
the evaluation of the derivative d(3{r-n., R)/dR: the most obvious choice for /5(r, R) would 
be a step- like function (tu = 0), but then eq. (9) becomes —2 A/? 5{rn — R)S{n). On 
a 3+1 dimensional lattice this clearly puts different values of i? to a completely unequal 
footing, and the whole operator makes sense only when one integrates it over some 
range of R. In principle this requires separate simulations for each value of R = in 
this range, which would be exceedingly costly. The form of fi{r, R) in eq. (8) smears 
the step over a shell R — w<rn^R + w. The larger w is, the smoother one can expect 
the behavior of F as a function of R to be, and F'{R) measured at any single value of 
R becomes meaningful. However, too large a w does not pin down the physical surface 
well enough. 

We chose w = 0.5a in order to keep the cross-over in /3(r, R) within one lattice spacing, 
which is anyway the shortest length scale on the lattice. It turns out that this choice 
makes the function F'{R) sufficiently regular (although there still remains some sys- 
tematic lattice effects, as seen from the similarity of the fluctuations of the data around 
the fits in Fig. 2.). 

We stress that the use of the w = instead of non-zero w does not make the calculations 
any more physical: when A(3 — > 0, the explicit w-dependence of the action vanishes, and 
the configurations extrapolate to physical bubble configurations. One can also expect 
that the effects of a particular choice of w become minimal. As discussed above, the 
bubble radius can be defined in numerous ways, and the discontinuity (or the cross-over) 
radius of the 'external field' A/?(r, R) is not particularly good parametrisation, since it 
is not any intrinsic property of the bubble. For example, it will turn out that the radius 
where the order parameter reaches the value half-way between the values in the interior 
and exterior bulk states can deviate from the cross-over radius considerably more than 
the half- width of the cross-over itself (Fig. 1). This implies that even with finite AjS 
the effects caused by choosing w = or w = 0.5 are overshadowed by the difference 
between the 'intrinsic' radius and the 'external field' radius of the bubble. We remind 
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again that this problem is only alleviated by the large value of the curvature term when 
compared to the interface tension. 

On the lattice we actually measure R in units of the lattice spacing a = 1/ (NtT) , and 
have, from eq. (4), 

djF/T) Att p,{T) - pf,{T) 8na{T)R 1 ^{T) 

R/a Nf T4 ^ a Nt ^ ^ 

Denoting the LHS by F'{R) we fit the lattice data to (from now on R/a — > R) 

F'{R) = csR' + C2R + ci, (11) 

extrapolate the coefRcients Cj to A(3 — and obtain the physical quantities at by 
comparison with eq. (10). 

We have performed simulations on a 2x 16^ lattice using an algorithm for SU(3) gauge 
theory with a combination of pseudo-heat bath and overrelaxation sweeps. The values 
of P were An/out = 5.1 ± A/3 with A/3 = 0.5, 0.25, 0.125 and the radius R varied between 
2 and 7 in steps of 0.5 with some intermediate points in steps of 0.25. The total number 
of iterations thus was 2(for H,Q)x 47(for values of A/3 and i?)x 10000 (or sometimes 
5000) iterations per point, amounting to a total of about IM iterations or (with 2.7 
s/sweep) about 700h of Cray XMP CPU time. By modern standard, the number of 
iterations on each point is relatively low, but the systems are rather far from the critical 
coupling and the autocorrelation times are small, <10 sweeps. 

The size of the lattice limits the physical bubble radius to < 3.5 /T. This is a serious 
limitation; it forces us to use relatively large A/3 to maintain the bubble on the lattice 
and to keep the correlation length small enough. Using eq. (6), one can see that at 
T = Tc the expected stable radius is 4/Tc and the radius where F{R) = F{0) = 
is 9/Tc. To truly study this region, the spatial size of the lattice should be at least 
40a. The small size of the bubbles and large A/3 values can cause systematical errors, 
whose magnitude remains unknown; for example, the terms 0{R^), n < in eq. (4) can 
become significant. Also, for a more detailed analysis of the A/3 dependence one would 
need more than 3 values of A/3 for both Q and H bubbles. 

Fig. 1 shows how the order parameter varies as a function of the distance from the 
center of the bubble for A/3 = 0.25; the pattern is similar for A/3 = 0.5 or 0.125. One 
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sees clearly the formation of an H bubble but with an interface thickness of some 4 units 
in R. In physical units the interface would be about 2 fm thick: we have 1 /Tg = N^a ^ 
1 fm, so that a ^ 0.5 fm. Remember that the input thickness w was only half a lattice 
spacing. 

The results for F'{R) and for the fitted parameters Ci are shown in Figs. 2 and 3 for H 
bubbles in Q and in Figs. 4 and 5 for Q bubbles in H. 

Since the numerical values of c^R^ are rather large, we shown in Fig. 2 the measured 
values of F'{R) after subtracting the term c^R^ with the fitted values of C3 shown in 
Fig. 3(a). Only values i? > 3 are used in the fit. The linear behaviour of the data is 
evident, although the x^/d.o.f. is as large as 2-4 depending on the lower cutoff in R. 

We attribute the large value of x^/d.o.f. to a lattice discretisation effect: the number 
of plaquettes inside the bubble does not grow smoothly with increasing R. 

The values of C3 in Fig. 3(a) extrapolate linearly to -0.38±0.10 when A/3 — > 0. Relative 
to the values extrapolated from, about -80,-40,-20, this is very close to the value 
expected at T = Tc (eq. (4)). Also the values for A/3 7^ are compatible with the 
pressure difference between homogeneous simulations at /3out and /3in. Similarly, C2 = 
(cq. (10)) extrapolates to the value 0.63±0.89 when A/3 — > 0. This coincides 
with the value (3) of a(Tc), but due to the large statistical errors we can confirm only 
the smallness of this value. Obviously, the methods utilised in [2,3] for planar interfaces 
are superior to the one used here for determining the interface tension a. The values of 
C2 for A/3 are discussed below. 

The main result is shown in Fig. 3(c): the curvature parameter ci = ^/{NfT"^) as given 
by the intercept of the lines in Fig. 3(a) at i? = 0. It is independent of A/3 and well 
compatible with the bag model result. The errors are of the same magnitude as with 
C2; however, due to the large value of ci the result is 3a away from zero. 

The corresponding data and fits for Q bubbles in H presented in Figs. 4 and 5 show a 
similar average pattern, but the error bars in the fits are so large that they make the 
result only la away from zero. In spite of this, the central value of ci coincides with 
that obtained for H bubbles, but has an opposite sign. Due to the large errors in Q 
bubbles, in what follows we concentrate mainly on the H bubble case. 
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From the derivative F'{R) we can obtain F{R) by integration. The result is shown in 
Fig. 6, normahsed arbitrarily to -F(O) = 0. One sees clearly the large numerical values 
for A/? ^ caused by the dominant term. The distance of extrapolation to A/? = 
is evidently large and we have constrained C3(A/? = 0) = 0. The resulting curve is in 
agreement with the phenomenological one in eq. (6). Notice, however, that one should 
not take the 1/R expansion too seriously for very small R. 

As discussed above, the definition of the physical radius of the bubble is necessarily 
ambiguous. It is essential to check the stability of the results against these reparametri- 
sations. Let us perform a reparametrisation of R consistent with the -R — > oo behaviour: 

R^ R + xi+X2/R+0{l/R^). (12) 

Then 

F{R) = \ czR^ + I (c2 + 2czXi)R^ + [ci + csa^i + cz{x\ + X2)]R + (13) 

We are mainly interested in the term linear in R. In it C3 is large but vanishes for 
A/3 and thus does not affect the physical value of ci. In the remaining term we 
know C2 from the value of the interface tension but cci, which corresponds to a constant 
shift in the radius, remains a free parameter. Even if we adjust it by ~ ±1 lattice units, 
which should encompass all 'reasonable' radius redefinitions discussed above, the linear 
term in R is only modified by ~ 10%. This is well below the statistical errors of ci. 

It is also possible to fix xi by comparing the second term with the values of the area 
term in planar interface simulations, i.e., the data in Fig. 3(b) with the data in Fig. 8 
of [2]. One obtains xi = 0.04 independent of A/3, which shows the consistency of the 
discussion. Thus at A/3 = the value of ci ^ —6 is only modified by a term of the 
magnitude 0.03; totally negligible. 

As discussed above, this seemingly surprising stability of ci is caused by the vanishing 
coefficient C3 when A/3 and by the very small value of the interface tension C2 when 
compared with ci. For A/3 the additional term 203X2 may explain why the data 
for ci seem independent of A/3 - clearly it is also possible that 7(T)/T^ is only weakly 
dependent on T, as the bag model predicts. 

In conclusion, our lattice Monte Carlo simulations support the existence of the curvature 
contribution 7 to the free energy of spherical bubbles predicted by the MIT bag model 
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[4]. The value predicted by the bag model is = ±327r/9, which is the same order of 
magnitude than our results; although for the quark bubbles the result is not statistically 
significantly non-zero. However, one should emphasise that the lattice used has a small 
Nt = 2, that even the smallest value of A/? still is quite large, that the bubble radii 
accessible to us may be too small to justify the truncation of the 1/R -expansion used 
here, and that there is no control of finite V effects. It is quite conceivable that more 
refined analyses will lead to a different numerical value, but we nevertheless find it 
quite interesting that Monte Carlo computations give such a clear signal for the MIT 
bag model value for the curvature term. For the interface tension a we could not obtain 
a statistically significant result. With further Monte Carlo runs one should also be 
able to make statements concerning F{R) when i? — and to address the question of 
stability of H bubbles for T > T^. quantitatively. 
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Figure Captions 



Fig. 1. The behaviour of the order parameter, the Wilson line L, as a function of the 
distance r from the center of the H bubble for A/3 = 0.25 and for i? = 3, 4, 5, 6, 7. The 
dotted line shows P{r,R = 5). 

Fig. 2. F'{R) - czB? for H bubbles. From up to down, A/3 = 0.5,0.25,0.125; the 
straight lines are the fits C2R + ci. 

Fig. 3. The fitted values of (a) C3, (b) C2, (c) function of A/3 for H bubbles with 

a linear fit to the data points and an extrapolation to A/3 = 0. In (c) the arrow shows 
the bag model point -167r/9. 

Fig. 4. F'{R) - csR'^ for Q bubbles. From up to down, A/3 = 0.5,0.25,0.125; the 
straight lines are the fits C2R + ci. 

Fig. 5. The fitted values of (a) C3, (b) C2, (c) function of A/3 for Q bubbles with 

a linear fit to the data points and an extrapolation to A/3 = 0. In (c) the arrow shows 
the bag model point lQn/9. 

Fig. 6. The free energy F{R)/T at T = Tc, normalised to F(0) = 0, integrated from 
the measured values of F'{R). R is in units of a and the curves are from below for 
A/3 = 0.5, 0.25, 0.125 and for the extrapolation to A/3 = 0. 



